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ABSTRACT 

We investigate the dynamical instability of the one-armed spiral m — I mode in differentially rotating 
stars by means of 3 -t- 1 hydrodynamical simulations in Newtonian gravitation. We find that both a 
soft equation of state and a high degree of differential rotation in the equilibrium star are necessary to 
excite a dynamical m = 1 mode as the dominant instability at small values of the ratio of rotational 
kinetic to gravitational potential energy, T/|W^|. We find that this spiral mode propagates outward from 
its point of origin near the maximum density at the center to the surface over several central orbital 
periods. An unstable m = 1 mode triggers a secondary 111 — 2 bar mode of smaller amplitude, and 
the bar mode can excite gravitational waves. As the spiral mode propagates to the surface it weakens, 
simultaneously damping the emitted gravitational wave signal. This behavior is in contrast to waves 
triggered by a dynamical m = 2 bar instability, which persist for many rotation periods and decay only 
after a radiation-reaction damping timescale. 

Subject headings: Gravitation — hydrodynamics — instabilities — stars: neutron — stars: rotation 
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1. INTRODUCTION 

Stars in nature are usually rotating and may be sub- 
ject to nonaxisymmetric rotational instabilities. An exact 
treatment of these instabilities exists only for incompress- 
ible equilibrium fluids in Newtonian gravity, (e.g. Chan- 
drasekhar 1969; Tassoul 1978). For these configurations, 
global rotational instabilities may arise from non-radial 
toroidal modes e™"^ (where m = ±1, ±2, . . . and ip is the 
azimuthal angle). 

For sufficiently rapid rotation, the m = 2 bar mode 
becomes either secularly or dynamically unstable. The on- 
set of instability can typically be identified with a critical 
value of the non-dimensional parameter /? = T/|VF|, where 
T is the rotational kinetic energy and W the gravitational 
potential energy. Uniformly rotating, incompressible stars 
in Newtonian theory are secularly unstable to bar-mode 
formation when f3 > Psec — 0.14. This instability can 
grow only in the presence of some dissipative mechanism, 
like viscosity or gravitational radiation, and the associ- 
ated growth timescale is the dissipative timescale, which 
is usually much longer than the dynamical timescale of the 
system. By contrast, a dynamical instability to bar- mode 
formation sets in when /3 > /3dyn — 0.27. This instabil- 
ity is independent of any dissipative mechanisms, and the 
growth time is the hydrodynamic timescale. 

Determining the onset of the dynamical bar-mode in- 
stability, as well as the subsequent evolution of an unsta- 



ble star, requires a fully nonlinear hydrodynamic simu- 
lation. Simulations performed in Newtonian gravity, (e.g. 
Tohline, Durisen, & McCoUough 1985; Durisen et al. 1986; 
WiUiams & Tohhne 1988; Houser, Centrella, & Smith 
1994; Smith, Houser, & Centrella 1995; Houser & Cen- 
trella 1996; Pickett, Durisen, & Davis 1996; Toman et al. 
1998; New, Centrella, & Tohhne 2000) have shown that 
Pdyn depends only very weakly on the stiffness of the equa- 
tion of state. Once a bar has developed, the formation of 
a two-arm spiral plays an important role in redistribut- 
ing the angular momentum and forming a core-halo struc- 
ture. Both /3dyn and /Sscc are smaller for stars with high de- 
gree of differential rotation (Tohline & Hachisu 1990; Pick- 
ett, Durisen, & Davis 1996; Shibata, Karino, & Eriguchi 
2002, 2003). Simulations in relativistic gravitation (Shi- 
bata, Baumgarte, & Shapiro 2000; Saijo et al. 2001) have 
shown that (3dyn decreases with the compaction of the star, 
indicating that relativistic gravitation enhances the bar 
mode instability. In order to efficiently use computational 
resources, most of these simulations have been performed 
under certain symmetry assumptions {e.g. 7r-symmetry) , 
which do not affect the growth of the m = 2 bar mode, 
but which suppress any m = 1 modes. 

Recently, Centrella et al. (2001) reported that such 
m — 1 "one-armed spiral" modes are dynamically un- 
stable at surprisingly small values of T/|iy|. Centrella 
et al. (2001) found this instability in evolutions of highly 
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differentially rotating equilibrium polytropes with poly- 
tropic index n = 3.33. Typically, these equilibria have a 
"toroidal" structure, so that the maximum density is not 
located at the geometric center but rather on a toroid ro- 
tating about the center. 

It is possible that the m = 1 instability in equilibrium 
stars is related to that arising in protostellar disk systems. 
This instability originally was found in nearly Keplerian, 
thin, gaseous disks around central point masses, both nu- 
merically (Adams, Rudcn, & Shu 1989; Hccmskcrk, Pa- 
palozizou, & Savonije 1992) and analytically (Shu et al. 
1990). The central point mass moves away from the cen- 
ter of mass of the whole system due to a perturbation 
and this displacement triggers the instability. This par- 
ticular mode of instability only occurs when the mass ra- 
tio Mdisk/-^totai exceeds 0.2. An m = 1 instability has 
also been found in thick, self-gravitating, protostellar tori 
(Woodward, Tohline, & Hachisu 1994) and protostellar 
disks (Laughlin & Bodenheimer 1994), as well as in finite 
fluid cores surrounded by disk halos (Pickett, Durisen, & 
Davis 1996). In the latter case, the instability arises from 
the internal interaction between different regions of a sin- 
gle, continuous body, and the disk does not need to satisfy 
the above mass criterion to trigger the m = 1 instabil- 
ity (see also Bonnell 1994, for an example of an unstable 
central accreting object surrounded by a rotationally sup- 
ported gas disk). 

The purpose of this paper is to study further the condi- 
tions under which a dynamical m = 1 instability is excited. 
We vary both the polytropic index, i.e. the stiffness of the 
equation of state, and the degree of differential rotation to 
isolate their effects on the instability. Since the onset of 
rotational instabilities is often characterized by /? we keep 
this value approximately fixed in our sequences. We find 
that a soft equation of state and a high degree of differen- 
tial rotation arc both necessary to dynamically excite the 
m = 1 mode at the small value of /3 = 0.14 chosen in this 
paper. We find that a toroidal structure is not sufficient 
to trigger the m = 1 instability, but our findings suggest 
that a toroidal structure may be necessary. 

While our goal is to gain a deeper understanding of the 
nature of the m = 1 instability as opposed to simulating 
realistic astrophysical scenarios, we point out that there 
exist evolutionary sequences that may well lead to rapidly 
and highly differentially rotating configurations. For ex- 
ample, cooling by thermal emission from a rotating star 
will cause the star to contract and spin up. If internal 
viscosity and magnetic fields are sufficiently weak, this 
process will lead to differential rotation even if the ini- 
tial configuration is rotating uniformly. This scenario may 
arise in supermassive stars, where the equation of state is 
dominated by radiation pressure and may be modeled by 
a (soft) n = 3 polytrope. In the absence of viscosity and 
magnetic braking, the star will contract quasi-statically 
as it cools to a toroidal configuration, which may be sub- 
ject to m = 1 or m = 2 instabilities (New & Shapiro 
2001). Stellar collisions and mergers may also lead to dif- 
ferentially rotating stars. For the coalescence of binary 
neutron stars (Shibata & Uryii 2000, 2002), the presence 
of differential rotation may temporarily stabilize the "hy- 
permassive" remnant and may therefore have important 



dynamical effects (Baumgarte, Shapiro & Shibata 2000; 
Lyford, Baumgarte & Shapiro 2003). However, as we find 
in this paper, the m = 1 mode is unstable only for very soft 
equations of state, so that it is not obvious that they will 
arise in the remnant of a binary neutron star merger. How- 
ever, they may arise in a rapidly spinning proto-neutron 
star core when surrounded by a fall-back disk, possibly 
forming a low mass condensation which can explode and 
induce a large neutron star recoil speed (Colpi & Wasser- 
man 2002). Finally, the m = 1 instability might arise in 
massive disks around black holes, especially if the disks are 
radiation-dominated, and hence governed by a soft equa- 
tion of state. 

This paper is organized as follows. In § 2 we present 
the basic equations, our initial data and diagnostics. We 
discuss our numerical results in § 3, and briefiy summa- 
rize our findings in § 4. Throughout this paper we use 
gravitational units ^ with G = c = 1 and adopt Cartesian 
coordinates {x,y,z). 

2. BASIC EQUATIONS 

2.1. Newtonian Hydrodynamics 

We construct a 3 -|- 1 dimensional Newtonian hydrody- 
namics code assuming an adiabatic F-law equation of state 
P = (F - l)pe, (1) 
where P is the pressure, F the adiabatic index, p the mass 
density and e the specific internal energy density. For per- 
fect fiuids the Newtonian equations of hydrodynamics then 
consist of the continuity equation 
dp d{pv^) 
dt dx^ 

the energy equation 



0, 



(2) 



.le-(r-i)p. ^ 
F 



(3) 



9^ 



(6) 



de d{ev^) 

'dt ^ dxi ~ 
and the Euler equation 

d{pvi) ^ dipvjV^) ^ 

dt dx^ 

Here is the fiuid velocity, $ is the gravitational poten- 
tial, and e is defined according to 

e = (p£)i/r. (5) 

We compute the artificial viscosity pressure -fvis from 
(Richtmyer & Morton 1994) 

CvisP{Sv)'^, for 5v < 0; 
0, for Sv > 0, 

where Sv = 2Sxdiv'^, Sx{— Ax ~ Ay = Az) is the local 
grid spacing and where we choose the dimensionless pa- 
rameter Cvis = 2. When evolving the above equations we 
limit the stepsize At by an appropriately chosen Courant 
condition. 

We have tested the ability of our code to resolve shocks 
by performing a wall-shock problem, in which two phases 
of a fluid collide at supersonic speeds. In Fig. 1 wc com- 
pare numerical results with the analytic solution for initial 
velocities that are similar to those found in our simulations 
below. With Cvis = 2 we find good agreement for Mach 
numbers up to Mmach < 6. The drop in density at a; = is 
usually interpreted as "wall heating" (e.g. Hawley, Smarr, 
& Wilson 1984). 



* Since we adopt Newtonian gravity in this paper, the speed of light only enters in the gravitational waveforms ( 



2.3 and § 3). 
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Fig. 1. — Comparison between numerical and analytical results in 
a one-dimensional wall shock problem at t = 2.5a;(''°""^'^/yo (where 
the fluid flow is aligned with the ii-axis). Solid and dashed Unes 
represent analytic and numerical results, respectively. For this sim- 
ulation we chose P = 1.30, = 2.80 x 10"^, k = 5.85 x lO'^, 
^(bound) ^ 0.5 with a grid space i5a; = 5 x 10-^ and vq = 2.78vs, 
where Vs is the initial speed of sound. 

The gravitational potential is determined by the Poisson 
equation 

A$ = Avrp. (7) 
which we solve subject to the outer boundary condition 



Here M is the total mass 
M = 



pdx^ 



and di is the dipole moment 

di = pxidx^. 
Jv 



(8) 
(9) 

(10) 



2.2. Initial Data 



As initial data, we construct differentially rotating 
equilibrium models with an algorithm based on Hachisu 
(1986). Individual models are parameterized by the ratio 
of the polar to equatorial radius i?p/i?oq, and a param- 
eter of dimension length d that determines the degree of 
differential rotation through 

Here $1 is the angular velocity, jo is a constant parameter 
with units of specific angular momentum, and w is the 
cylindrical radius. The parameter d determines the length 
scale over which changes; uniform rotation is achieved 
in the limit — > oo. For the construction of initial data 
we also assume a polytropic equation of state 

P = Kpi+i/", (12) 

where n=l/(r— l)is the polytropic index and k a con- 
stant. In absence of shocks, the polytropic form of the 
equation of state is conserved by the F-law equation of 
state (eq. (1)). 

To enhance any m = 1 or m = 2 instability, we disturb 
the initial equilibrium density Peq by a non-axisymmetric 
perturbation according to^ 



y 



R. 



eq 



-"-eq 



(13) 



with 5^^^ = 5^"^^ =10 ^ in all our simulations. 



2.3. Gravitational Waveforms 

We compute approximate gravitational waveforms by 
evaluating the quadrupole formula. In the radiation zone, 
gravitational waves can be described by a transverse- 
traceless, perturbed metric hf^^ with respect to a flat 
spacetime. In the quadrupole formula, /i^'^ is found from 
(Misner, Thorne, & Wheeler 1973) 



h 



TT 



2 d 

rdt^ ' 



(14) 



where r is the distance to the source, lij the quadrupole 
moment of the mass distribution (see eq. [36.42b] in Mis- 
ner, Thorne, & Wheeler 1973), and where TT denotes the 
transverse-traceless projection. Choosing the direction of 
the wave propagation to be along the z axis, the two po- 
larization modes of gravitational waves can be determined 
from 



2 V^xx 



and hy 



'TT 



For observers along the ^;-axis, we thus have 

lull I ^ 



r/i+ 1 d ■ 



M 



Mdi"""- 



(15) 

(16) 
(17) 



The number of time derivatives 1^ that have to be taken 
can be reduced by using the continuity equation (2) 

i^■ = J [pv'x^ + px'v^)d^x, (18) 

in equations (16) and (17) (see Finn 1989). 

2.4. Numerical Implementation and Diagnostics 

Our code is based on the post-Newtonian hydrodynam- 
ics scheme of Shibata, Baumgarte, & Shapiro (1998) and 
Saijo et al. (2001), to which the reader is referred for a 
more detailed description, discussion and tests. We choose 
the axis of rotation to align with the z axis, and assume 
planar symmetry across the equator. The equations of 
hydrodynamics are then solved on a uniform grid of size 
169 X 169 X 85. We terminate our simulations either when 
the central density has increased to a point at which our 
resolution becomes inadequate, or after a sufficient number 
of central rotation periods (between 20 and 40) in order 
for us to detect dynamical instabilities. 

We monitor the conservation of mass M (eq. [9]), angu- 
lar momentum J 



J = y p{xvy - yv'^)d^x, 



(19) 



energy E 



E = T+U + W =^ J pv'vid^x + j ped^x+^j p^d^ 



(20) 



^ The numerical finite difference error is in principle sufficient to trigger instabilities, but starting from such a small amplitude it would taJje 
the instability prohibitively long to reach saturation. 
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Table 1 

Initial data for bar formation tests {n = 1). 



d/Rcq ^ 






Pc / Pmax ^maxd/-^cq 


T/\W\* m = l 


m = 2 


0.20 


0.250 


26.0 


0.160 0.383 


0.119 Stable 


Unstable 



° i?oq = equatorial radius. 
* i?p = polar radius. 

ric = central angular velocity, fioq = equatorial angular velocity at the surface. 

Pc — central density; pmax = maximum density. 

^maxd = distance between the origin and the location of maximum density. 
^ T = rotational kinetic energy; W — gravitational potential energy. 



and the location of the center of mass 



'CM 



^CM 



(21) 



Here T is the kinetic energy (all rotational a,t t = 0), U 
the internal energy, and W the gravitational potential en- 
ergy. Given our assumption of equatorial symmetry, we 
have X 



CM 



identically, so that we only need to mon- 
itor the X and y components of x^y^. Due to our flux- 
conserving difference scheme the mass M is also conserved 
up to round-off error, except if matter leaves the compu- 
tational grid (which was never more than 0.01 % of the 
total mass). In all cases reported in § 3 the energy E and 
the angular momentum J were conserved up to ~ 0.1% of 
their initial values, and the center of mass moved by less 
than about 1% of one spatial grid cell per central rotation 
period. 

To monitor the development of m = 1 and 171 — 2 modes 
we compute a "dipole diagnostic"^ 

D = (^™^)_r - I P^^FAd^- (22) 



and a "quadrupole diagnostic" 



^ (x^-,J) + .(2.,) ^3,^ (23) 



where a bracket denotes the density weighted average. In 
the following we only plot the real parts of D and Q. 

3. RESULTS 

3.1. Dynamical bar formation 

Before studying m = 1 one-armed spiral instabilities, it 
is useful to test the capability of our code and our diag- 
nostics to detect any instabilities. To do so, we reproduce 
an m = 2 bar mode instability that was recently found by 
Shibata, Karino, & Eriguchi (2002) in highly differentially 
rotating n — \ polytropes for surprisingly small values of 
T/|VF|. The parameters of our initial data for this test are 
listed in Table 1. For all our simulations we add a small 
dipole (m — 1) and quadrupole (m = 2) perturbation to 
the initial equilibrium star (eq. [13]) to enhance the growth 
of any instability. 
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Fig. 2. — Diagnostics D and Q as a function of t/Pc for our bar 
formation model (see Table 1). Solid and dotted lines denote the 
values of D and Q, respectively. We terminate our simulation when 
t = 132Pc. Hereafter Pc represents the central rotation period. 
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Fig. 3. — Final density contours in the equatorial plane for our 
bar formation model. Snapshots are plotted at (t/Pc, Pmax/pSax) 
= (132, 1.25), where 

pmax is the maximum density, Pmax is the ini- 
tial maximum density, and R is the initial equatorial radius. The 
contour lines denote densities p/pmax = 6.67 X (16 -i) X 10-2 (i = 
1,---,15). 



^ Our diagnostics differ from those in previous treatments (CentrcUa et al. 2001), where the growth of the mode is measured at a single, 
arbitrary Eulerian radius in the equatorial plane inside the star. 
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Fig. 4. — Gravitational waveform for bar-unstable star as seen 
by a distant observer located on the z-axis. 

In Fig. 2 wc show both diagnostics D and Q as a func- 
tion of time. The dipole diagnostic D remains very small 
throughout the evolution (small oscillations are due to the 
initial perturbation), while the quadnipole diagnostic Q 
grows exponentially until it saturates. These results indi- 
cate that the star is unstable towards bar formation but 
stable towards one-armed spiral formation. The bar per- 
sists without decay for over one surface-rotation period 
following saturation, corresponding to over 30 central ro- 
tation periods. After this we terminate our integration. 

The bar mode formation is also evident in Fig. 3, which 
shows a snapshot of the density contours just before we ter- 
minate the evolution. Owing to the small value of T/|W| 
the bar is too weak to form double spiral arms. The grav- 
itational waveform emitted by the bar formation is shown 
in Fig. 4. We expect that it will survive without decay un- 
til gravitational radiation-reaction forces destroy the bar 

These simulations indicate that our code and diagnos- 
tics arc capable of detecting instabilities, and also recon- 



firm the findings of Shibata, Karino, & Eriguchi (2002) 
that strongly differentially rotating stars can be unstable 
to dynamical bar mode formation even at very small values 
oiT/\W\. 



3.2. Dynamical one-armed spiral formation 

We now focus on m = 1 one-armed spiral instabilities. 
Before we analyze their dependence on the stiffness of the 
equation of state and the degree of differential rotation in 
the following subsections, we first want to reconfirm the 
findings of Centrella et al. (2001). To reconstruct their 
initial data, we adopt a polytropic index of n = 3.33 and a 
high degree of differential rotation {d/Re = 0.2). We study 
two different models, which are detailed in Table 2. The 
more rapidly rotating Model I (a) (the case r/|I^| = 0.14 
of Centrella et al. 2001) has a toroidal structure, while 
Model I (b) (the case T/\W\ = 0.09 of Centrella et al. 
2001) does not. Confirming the results of Centrella et al. 
(2001), we find that Model I (a) develops an m = 1 insta- 
bility, while Model I (b) remains stable. 

The different stability properties of the two models can 
be seen in Fig. 5, where we show both diagnostics D and 
Q. For Model I (b), both diagnostics remain very small, 
indicating stability', while for Model I (a) both diagnos- 
tics grow. The dipole diagnostic D, however, grows more 
strongly than the quadrupole diagnostic Q, indicating that 
the m = 1 mode is the dominant unstable mode. This is 
also evident in the density contours in Fig. 6, which clearly 
exhibit the one-armed spiral in Model I (a) at intermedi- 
ate times. In all cases that we found to be unstable to an 
TO = 1 mode, we simultaneously found a growing to = 2 
mode. 

In Fig. 7 we show the maximum density Pmax as a func- 
tion of time for both models. Even for the stable Model I 
(b) the central density slowly increases over the course of 
several central rotation periods. This slow growth is due to 



^ The small growth of the m = 1 mode in Model I (b) is a numerical artifact triggered by the initial perturbation; the absence of an exponential 
growth indicates that this is not an instability. 
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Fig. 5. — Diagnostics D and Q as a function of t/Pc for Model I (a) and (b) (see Table 2). Solid and dotted lines denote D and Q. We 
terminate our simulation at t ~ 20Pc or when the maximum density of the star exceeds about 10 times its initial value. 
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Table 2 

Initial data for the m — 1 test (n = 3.33). 



Model 


d/ Rcq 


-^p/-^cq 


/ ^cq 


Pc/Pmax 




T/\W\ 


m = 1 


1(a) 


0.20 


0.417 


26.0 


0.531 


0.192 


0.144 


Unstable 


1(b) 


0.20 


0.542 


26.0 


1.00 


0.00 


0.090 


Stable 





- 1.5- 1 - 0.5 0.5 1 1.5 - 1.5- 1 - 0.5 0.5 1 1.5 

X / R X / R 

Fig. 6. — Intermediate and final density contours in the equatorial plane for Model I (a) and Model I (b). Snapshots are plotted at (t/Pc, 
Pmax/PmL, d) = (a)-i (16.3, 3.63, 0.287), and (b)-i (14.7, 2.08, 0.333), (a)-u (23.3, 11.5, 0.287), and (b)-ii (20.6, 3.66, 0.333). The contour 
lines denote densities p/pmax = lO"'^^"''"^ (i = 1, • ■ ■ , 15). 



numerical and artificial viscosity, which tends to decrease 
the degree of differential rotation. As a consequence, the 
angular velocity at the center decreases, which also de- 
creases the rotational support of the matter at the center, 
and hence leads to a slow increase of the central density, 
even for supposedly stable stars (see also Fig. 8). This ef- 
fect is a numerical artifact, although viscosities in stars in 
nature would have a very similar effect. For the unstable 
Model I (b), however, we find a much more rapid increase 
in the central density. This enhanced increase is caused by 
the growing spiral instability as it redistributes the matter 
in the star and destroys the toroidal structure (compare 
Fig. 8). 

Unlike for bar formation, where the bar persists for 
many rotational periods (compare § 3.1, Brown 2000; Saijo 



et al. 2001), we find that D and Q start decreasing imme- 
diately after reaching a maximum (see Fig. 5, note that 
the decrease in D is not as dramatic as the decrease in 
Q). This is also evident in Fig. 6, where the density con- 
tours approach axisymmetry at late times. As the spiral 
arm propagates through the star, it rearranges the density 
profile, eliminates the toroidal structure, and ultimately 
leads to a new axisymmetric equilibrium configuration. 
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Fig. 7.- Maximum density pmax as a function of t/Pc for Model 
I (a) (solid) aud JModcl 1 (b) (dotted). We terminate our simulation 
at t ~ 20Pc or when the maximum density of the star exceeds about 



10 times its initial value p, 



(0) 




Fig. 8. — Density profiles along the x-axis during the evolution for 
Models I (a) and I (b). Solid, dotted, dashed, dash-dotted lines de- 
note times i/fc = (a) (1.16x10-^, 6.99, 14.0, 21.0), (b) (7.36x10-", 
6.63, 13.3, 19.9), respectively. Note that the density distribution de- 
velops asymmetrically in the presence of the m = 1 mode instability, 
and that this instability destroys the toroidal structure. 
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Flc;. 9. — Gravitational waveforms as seen by a distant observer 
located on the 2;-axis for Model I (a) (solid line) and Model I (b) 
(dashed line). 

In Fig. 9 we show the gravitational wave signal emitted 
from this instability. Gravitational radiation couples to 
quadrupole moments, and the emitted radiation therefore 
scales with the quadrupole diagnostic Q, which wc always 
find excited along with the m = 1 instability. We con- 
sistently find that the pattern period of the the m = 2 
modes is very similar to that of the m = 1 mode, suggest- 
ing that the former is a harmonic of the latter (see Table 
3). Since the diagnostic Q does not remain at its maximum 
amplitude after saturating, we find that the gravitational 
wave amplitude is not nearly as persistent as for the bar 
mode instability. We also find that the gravitational wave 
period, here Pgw ^ 0.7Pc ~ i is different from the 
value Pgw ~ 3.3Pc ~ we found for the bar mode in 
§ 3.1, which points to a difference in the generation mecha- 
nism. Characteristic wave frequencies faw correspond to 
the central rotation period of the star. 
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t/P t/P 



Fig. 10. — Diagnostics D and Q as a function of t/Pc ^or Models II (see Table 4). Solid and dotted lines denote D and Q. We terminate 
our simulation at t ~ 25Pc or when the maximum density of the star exceeds about 10 times its initial value. 
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- 1.5- 1 - 0.5 0.5 1 1.5 - 1.5- 1 - 0.5 0.5 1 1.5 

X / R X / R 



Fig. 11. — Intermediate and final density contours in the equatorial plane for Models II. Snapshots are plotted at (t/Pc, Pmax/Pmaxi d-) = 
(a)-i (23.8, 1.14, 0.220), (b)-i (20.6, 1.90, 0.220), (c)-i (17.3, 4.98, 0.287), (d)-i (16.3, 3.63, 0.287), (a)-ii (34.7, 1.24, 0.220), (b)-ii (30.1, 2.92, 
0.220), (c)-ii (25.2, 7.41, 0.287), and (d)-ii (23.3, 11.5, 0.287). The contour lines denote densities p/pmax = 10'^'^'^-^^'^ {i = 1, ■ • • , 15). 
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Table 3 

Comparison of Pattern Periods 



II (b) b 
II (c) b 

II (d) b 

III (c) 
III (d) 



Model Period of m = 1 mode Pc Period of m = 2 mode Pc Pattern Period Pr 



0.7 
1.2 
1.0 
0.7 
0.7 
0.7 



1.6 
2.5 
2.0 
1.6 
1.6 
1.6 



0.8 
1.3 
1.0 
0.8 
0.8 
0.8 



" Sec Fig. 5 
^ See Fig. 10 
^ See Fig. 15 




t/P 



The results of this Subsection confirm the findings of 
Centrella et al. (2001), and establish that stars with soft 
equations of state and large degrees of differential rota- 
tion are unstable to one-armed spiral arm formation. Such 
stars have a toroidal structure which is erased by the grow- 
ing m = 1 mode. One might be lead to believe that this 
toroidal structure is a necessary and perhaps even a suf- 
ficient condition for the growth of the m = 1 instability. 
In the following two Subsections we analyze the depen- 
dence of the onset of instability on both the stiffness of 
the equation of state and the degree of differential rota- 
tion, and find that toroidal structure alone is not sufficient 
for a one-armed spiral instability. 



Fig. 12. — Maximum density pmax as a function of t/Pc for Model 
II (a) (solid line), Model II (b) (dotted line), Model II (c) (dashed 
line), and Model II (d) (dash-dotted line). 



3.3. Stiffness of the equation of state 

We parameterize the stiffness of the equation of state 

by varying the polytropic index n between n = 3.33 and 
n = 2. In this sequence we keep the degree of differential 




x/R x/R 

Fig. 13. — Density profiles along the x-axis during the evolution of Models II. Solid, dotted, dashed, dash-dotted line denote at time t/Pc = 
(a) (2.17 X 10-3, 10.8, 21.7, 32.5), (b) (1.56 x 10"^, 9.56, 19.0, 28.5), (c) (1.33 x 10"^, 7.97, 15.9, 23.9), (d) (1.16 x 10"^, 6.98, 14.0, 21.0), 
respectively. Note that the toroidal structure vanishes at late times for Models II (b), II (c), and II (d). 
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Fig. 14. — Gravitational waveforms as seen by a distant observer located on the ^-axis for Models II. 
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Fig. 15. 
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t/P t/P 

c c 

Diagnostics D and Q as a function of t/Pc for Models III (see Table 5). Solid and dotted lines denote D and Q. 



ONE-ARMED SPIRAL INSTABILITY IN DIFFERENTIALLY ROTATING STARS 



11 



pi, 



PS 




1.5 - 1 - 0.5 0.5 1 1.5 

X / R 



1.5 - 1 - 0.5 0.5 1 1.5 

X / R 



Fig. 16. — Intermediate and final density contours in the equatorial plane for Models III. Snapshots are plotted at (t/Pc, Pmax/pmaxi <i) = 
(a)-i (6.81, 1.01, 0.200), {b)-i (10.4, 1.07, 0.267), (c)-i (17.6, 3.72, 0.287), (d)-i (16.3, 3.63, 0.287), (a)-ii (9.90, 1.02, 0.200), (b)-ii (14.8, 1.09, 
0.267), (c)-ii (25.6, 5.05, 0.287), and (d)-ii (23.3, 11.5, 0.287). The contour lines denote densities p = pmax x 10-(lS-')'*(i = 1, ■ ■ ■ , 15). 
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Table 4 

Initial data sequence varying the polytropic index. 



Model 




d/Rcq 


-^p/-^cq 


/ ^cq 


Pc/ Pmax 


-^maxd / -^cq 


T/\W\ 


m = 1 stability 


II (a) 


2.00 


0.20 


0.271 


26.0 


0.091 


0.349 


0.145 


Stable 


II (b) 


2.50 


0.20 


0.354 


26.0 


0.193 


0.295 


0.145 


Unstable 


II (c) 


3.00 


0.20 


0.396 


26.0 


0.325 


0.243 


0.147 


Unstable 


II (d)b 


3.33 


0.20 


0.417 


26.0 


0.531 


0.192 


0.144 


Unstable 



° polytropic index 
'' Same as Model I (a). 



rotation (i.e. d and hence Qc/^cq) fixed, and adjust the 
overall rotation rate (parameterized by Rp/Rcq) so that 
the value of remains very close to 0.144 (as for 

Model I [a]). We list our four different Models II in Table 
4, and note that Model II (d) is identical to Model I (a). 

Figure 10, where we plot the dipole diagnostic I? as a 
function of time, clearly shows that an m = 1 instability is 
excited in Models II (b) and II (c) in addition to Model II 
(d). After reaching saturation, D decreases again, similar 
to Model I (a) which we described in detail in § 3.2. Model 
II (a), however, which has the most pronoimced toroidal 
structure, remains stable. These findings are also evident 
in Fig. 11, where we show density contours of intermediate 
and final configurations. 

In Fig. 12, we show the maximum density as a func- 
tion of time. As we have seen in § 3.2, the maximum 
density slowly increases in all cases due to dissipation of 
differential rotation. Once the one-armed spiral forms in 
Models II (b) and II (c), however, this increase is much 
more rapid, which indicates again that the unstable mode 
rearranges the matter in the star and destroys the toroidal 
structure. This effect can also be seen in the density pro- 
files in Fig. 13. 

We show the gravitational wave signal emitted from 
Models II in Fig. 14. As we found in § 3.2, and consis- 
tent with the diagnostics D and Q, the gravitational wave 
signal emitted by the one-armed spiral mode does not per- 
sist over many rotational periods, and instead decays fairly 
rapidly after it has been excited. This characteristic is very 
different from what has been found for m = 2 bar mode 
instabilities (compare § 3.1 and Brown 2000; Saijo et al. 
2001). We also find that the maximum wave amplitude is 
much smaller than can be found for configurations unsta- 
ble to a pure bar mode (compare Fig. 4) as gravitational 
radiation requires a quadrupole distortion and the m = 2 
perturbation in Models II is only being excited as a lower- 



amplitude harmonic of the m = 1 mode (as suggested by 
a comparison of the pattern periods; see Table 3). 



10 



Q. 



10, 



10 



t/p 



I I I I L 

15 20 



25 



Fig. 17. — Maximum density pmax as a function of t/Pc for 
Model III (a) (solid line), Model III (b) (dotted line), Model III (c) 
(dashed line), and Model III (d) (dash-dotted line). 

3.4. Degree of differential rotation 

We now focus on the dependence of the one-armed spiral 
instability on the degree of differential rotation. Starting 
again with Model I (a), we now increase the parameter d to 
explore more modest degrees of differential rotations. As 
before, we would like to keep /? ^ 0.14 in this sequence. 
For very soft equations of state, this value can only be 
achieved for very strong degrees of differential rotation. 
Therefore, in order to keep (3 approximately constant, we 
simultaneously have to decrease n as we decrease the de- 



Table 5 

Initial data sequence varying the degree of differential rotation. 



d/Req Rp/Req i^c/^eq Pc/ P: 



/Reg T/\W\ m = l stability 



Model 
III (a) 

III (b) 

III (c) 

ni (.1) 



R, 



■maxd/ 



1.00 


0.62 


0.500 


3.60 


0.992 


0.189 


0.150 


Stable 


2.00 


0.41 


0.479 


6.95 


0.935 


0.198 


0.150 


Stable 


3.00 


0.25 


0.438 


17.0 


0.695 


0.197 


0.147 


Unstable 


^' :-!.:-!:-! 


0.20 


0.417 


2().() 


().5:-il 


0.192 


0.144 


Uust.ahlo 



Model I (a) in Table. 2. 
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Fig. 18. Density profiles along the x-axis during the evolution for Models III. Solid, dotted, dashed, dash-dotted line denote at time 
tjP^ = (a) (3.09 X 10--*, 3.09, 6.19, 9.28), (b) (4.95 x lO"", 4.95, 9.89, 14.8), (c) (8.82 x lO"*, 8.82, 17.6, 26.4), (d) (1.16 x IQ-^, 6.99, 14.0, 
21.0), respectively. Note that the toroidal structure vanishes at the late time in Models III (c) and III (d). 
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Fig. 19. — Gravitational waveforms as seen by a distant observer located on the z-axis for Models III. 
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gree of differential rotation*. We list the details of our 
Models III in Table 5. 

We show the dipole diagnostic D as a function of time 
in Fig. 15, which shows that Models III (a) and III (b) 
are stable against one-armed spiral formation while Mod- 
els III (c) and III (d) (which is the same as Model I [a]) 
are not. The same conclusion can be drawn from the den- 
sity snapshots in Fig. 16. As in § 3.3, we find that the 
one-armed spiral results in a large increase in the central 
density (Fig. 17), and an elimination of the toroidal struc- 
ture (Fig. 18). Tohline & Hachisu (1990) similarly found 
that the elimination of the toroidal structure is related to 
an outward transport of angular momentum. To quantify 
this effect, we monitored the angular momentum distri- 
bution in Model Ill(d) by computing a mean radius of 
angular momentum 

^ , ^ M[J dvpy^x^ + y'^{xvy - yv'^)] ^^^^ 
J[J dvp\J + y^] 

Initially this mean radius is 1.1, but increases to 1.5 at 
t = 24.5Pc, indicating that in fact the m = 1 mode trans- 
ports angular momentum outward. 

We show gravitational waveforms from Models III in 
Fig. 19. We again find that the amplitude decreases after 
reaching a maximum, which is a typical behavior of m = 1 
instability (see Table 3 that the pattern period of diag- 
nostics D and Q are the same). In some cases, however, 
this decrease is not monotonic, and the amplitude may in- 
crease again to form several distinct wave packets. Our 
numerical data are not sufficient to determine the generic 
character of the gravitational waves emitted from m = 1 
instabilities, and we expect that this will be subject of fu- 
ture investigations. The problem is that the growth of cen- 
tral concentration during the evolution exceeds the ability 
of our code to resolve the innermost regions for arbitrary 
long times in all cases. 

4. DISCUSSION 

We have studied the conditions under which Newtonian, 
differentially rotating stars are dynamically unstable to an 
m = 1 one-armed spiral instability, and found that both 
soft equations of state and a high degree of differential 
rotation arc necessary to trigger the instability. For suffi- 
ciently soft equations of state and sufficiently high degrees 



of differential rotation we found that stars are dynami- 
cally mistable even at the small values of T/IVI^I ~ 0.14 
considered in this paper. 

While we find that a toroidal structure alone is not suf- 
ficient for the m = 1 instability, all the models that are 
unstable do have a toroidal structure, suggesting that this 
may be a necessary condition. The growing m = 1 mode 
redistributes both matter and angular momentum in the 
unstable star and destroys the toroidal structure after a 
few central rotation periods. 

Quasi-periodic gravitational waves emitted by stars with 
m = 1 instabilities have smaller amplitudes than those 
emitted by stars unstable to the m ~ 2 bar mode. For 
TO = 1 modes, the gravitational radiation is emitted not 
by the primary mode itself, but by the m = 2 secondary 
harmonic which is simultaneously excited, albeit at a lower 
amplitude (see Fig. 5). Unlike the case for bar-unstable 
stars, the gravitational wave signal does not persist of 
many periods, but instead is damped fairly rapidly in most 
of the cases we have examined. 

We have plotted typical wave forms for stars unstable to 
TO = 2 bar modes in Fig. 4 and for stars unstable to one- 
armed spiral TO = 1 modes in Figs. 9, 14 and 19. Charac- 
teristic wave frequencies /gw are seen to be ~ ~ fici 
and are considerably higher than fieq ^ {M/R^Y^"^ due to 
appreciable differential rotation. For supermassive stars 
(M > 10^ Mq) the amplitudes and frequencies of these 
waves fall well within the detectable range of LISA (see, 
e.g., New & Shapiro (2001)). 
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